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ABSTRACT 

\j 

An  analytic  procedure  was  developed  to  predict  the  expected 
proportion  of  defectives  produced,  using  a  specified  set  of 
fabrication  and  assembly  processes.  This  procedure  was  demonstrated 
for  a  General  assembly  process,  a  Peg-in-Hole  assembly  process  and  a 
Box-and-Cover  assembly  process.  This  procedure  was  compared  to  the 
formulation  of  a  Monte  Carlo  simulation  model  and  shown  to  require 
the  same  initial  information.  Two  applications  of  the  results  of 
this  procedure  were  developed.  The  first  was  an  optimization  model 
used  to  allocate  a  fixed  amount  of  capital  improvements  to  a  discrete 
number  of  possible  candidate  assembly  processes.  The  second  was  to 
statistically  control  an  assembly  process  using  a  p  Chart.  A 
numerical  example  illustrating  the  procedure  and  its  applications  was 


also  presented. 
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Chapter  1 


INTRODUCTION 


Manufacturing  processes  are  generally  thought  of  as  the 
fabrication  of  component  parts  followed  by  assembly  into  a  final 
product.  Traditionally,  the  focus  of  Statistical  Process  Control  or 
SPC  has  been  to  control  the  fabrication  processes  of  the  component 
parts.  The  techniques  used  to  control  these  processes  (i.e.,  proces 
capability  studies,  control  charts,  etc.)  are  well  developed  and 
widely  used  in  industry.  Assembly  processes  are  also  stochastic  in 
nature  and  can  be  controlled  by  the  same  statistically  based 
techniques  that  are  used  for  fabrication  processes.  Since  the  final 
product  of  a  given  manufacturing  process  almost  always  involves 
assembly  operations,  the  ability  to  determine  if  the  assembly 
processes  are  in  statistical  control  would  be  very  helpful  in 
maintaining  the  quality  level  of  the  final  product.  In  addition, 
knowing  the  expected  number  of  defectives  for  a  given  set  of 
fabrication  and  assembly  processes  can  be  very  helpful  when 
allocating  resources  for  process  improvements. 

Several  literature  reviews  of  existing  research  in  the  areas  of 
assembly  process  variation  and  control  yield  four  basic  groups  of 
research.  The  first  group  of  papers  [6,7,12]  studies  the  problem  of 
assembly  process  variation  using  the  Monte  Carlo  simulation  method. 
While  this  technique  is  useful  to  get  a  general  idea  of  how  a 
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complete  manufacturing  process  operates,  it  does  not  lend  itself  to 
the  problem  of  optimally  allocating  capital  improvements  or  real  time 
process  control.  Another  group  of  papers  [8,9,10,11,14]  considers 
the  effect  of  assembly  errors  on  circular  arc  or  Novikov  types  of 
gears.  These  papers  consider  the  assembly  errors  to  be  deterministic 
and  measurable  for  a  given  gear  pair.  This  analysis  is  useful  for 
failure  mode  analysis  and  other  preliminary  design  related  tasks,  but 
not  for  true  manufacturing  and  quality  control  applications. 

The  final  two  groups  of  research  are  more  closely  related  to  the 
main  thrust  of  this  paper.  The  first  of  these  groups  is  research 
into  the  allocation  of  errors  in  hierarchical  calibration  or  assembly 
systems.  Rajaraman  [13]  expands  Crow's  work  [3,4]  on  calibration 
systems  to  assembly  systems.  Rajaraman  develops  a  dynamic 
programming  model  to  optimally  allocate  resources  to  each  stage  of 
assembly.  The  main  drawback  to  this  approach  is  that  the  assembly 
variations  at  each  stage  have  to  be  independent  to  use  a  dynamic 
programming  approach.  While  this  is  true  for  the  linear  or 
"telescopic"  assembly  process  described  in  [13],  it  isn't  true  for  a 
general  or  "tree-like"  assembly  process.  The  problem  given  in 
Chapter  3  of  this  paper  is  an  example  of  a  "tree-like"  assembly 
process  since  the  assembly  variations  in  the  third  stage  are 
dependent  on  the  variations  in  the  first  stage.  Rajaraman  also 
illustrates  that  the  manner  in  which  errors  propagate  between  stages 
is  unknown  and  can  be  assumed  to  be  a  linear  or  nonlinear  function. 

The  last  group  of  research  related  to  assembly  process  variation 
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is  the  analysis  of  basic  assembly  operations,  used  in  the 
implementation  of  robotic  or  automated  assembly  stations.  While 
there  are  a  large  number  of  papers  to  choose  from,  the  works  by 
Abdel-Malek  [1]  and  Boucher  [2]  are  used  for  their  clear  analyses  of 
the  two  most  basic  assemblies,  Peg-in-Hole  and  Box-and-Cover .  Both 
authors  develop  probabilistic  relationships  for  successful  completion 
of  both  of  these  assembly  operations.  Their  research  concentrates  on 
the  probability  of  a  robot  with  known  repeatability  being  able  to 
insert  a  peg  into  a  hole  or  a  box  into  a  cover.  They  do  not  consider 
the  geometries  of  the  actual  parts  that  are  assembled  and  whether 
these  parts  are  assembled  within  acceptable  tolerances. 

There  are  three  basic  objectives  of  this  paper.  The  first 
objective  is  to  analytically  predict  the  expected  number  of  defective 
units  for  a  given  set  of  fabrication  and  assembly  processes.  This 
will  be  accomplished  using  the  same  probabilistic  analysis  used  for  a 
Monte  Carlo  type  of  simulation  except  that  the  results  will  be 
determined  using  numerical  integration  and  not  simulation.  The  next 
objective  is  to  develop  an  optimization  model  to  allocate  a  limited 
amount  of  capital  resources  to  reduce  the  variation  in  the  existing 
assembly  processes.  Finally,  a  method  for  statistically  controlling 
an  assembly  process  is  developed. 


4 


Chapter  2 

ANALYTIC  METHODOLOGY 


The  basic  approach  used  in  this  method  is  to  analyze  a  given 
assembly  process  from  a  probabilistic  viewpoint.  There  are  three 
basic  probabilities  associated  with  any  assembly  process.  First, 
there  is  the  probability  that  assembly  is  possible.  For  example,  it 
is  possible  that  machined  components  will  not  fit  together  even  if 
there  is  no  assembly  process  variation.  The  next  event  of  interest 
is  the  probability  that  assembly  will  be  successful  in  light  of  the 
existing  assembly  process  variation.  These  two  probabilities  have 
been  most  widely  researched  for  application  in  robotic  assembly,  as 
mentioned  in  Chapter  1.  Finally,  there  is  the  probability  that  a 
successfully  assembled  unit  meets  the  specified  quality  requirements. 
To  have  a  completed  assembly  with  the  required  quality,  all  three  of 
these  probabilities  have  to  be  satisfied. 

As  stated  in  Chapter  1,  the  most  common  methodology  used  to 
obtain  an  expected  number  of  defective  assemblies  (i.e.,  an  assembly 
that  doesn't  meet  all  three  conditions  above)  is  the  Monte  Carlo 
Simulation.  In  this  paper,  an  analytic  solution  is  found  using 
numerical  integration  techniques.  The  method  used  to  find  this 
solution  will  be  demonstrated  on  three  basic  assembly  processes.  The 
first  is  a  General  assembly  process,  one  that  doesn't  contain  any 
assembly  aids  (i.e.,  fixtures,  pegs  and  holes,  etc.).  Since  the 
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assembly  variation  associated  with  this  process  isn't  limited  by 
assembly  aids,  the  assembly  variation  isn't  restricted  in  any 
direction.  Next,  a  Peg-in-Hole  assembly  will  be  considered.  The 
Peg-in-Hole  assembly  process  limits  the  possible  variation  in  the  X 
and  Y  directions.  Finally,  a  Box-and-Cover  assembly  will  be 
considered.  The  Box-and-Cover  assembly  process  is  essentially  a 
rectangular  Peg-in-Hole  which  limits  variation  in  the  X,  Y  and  9 
directions.  It  is  assumed  during  this  analysis  that  only  one  attempt 
is  made  at  assembly.  Also,  if  the  assembly  is  successful,  the  part 
is  fixed  to  the  work,  capturing  the  assembly  errors. 

Finally,  a  review  of  the  basics  of  numerical  integration 
including  the  definition  of  some  useful  notation  is  appropriate 
before  beginning  the  analysis.  The  general  process  is  to  write  the 
expression  of  interest  as  a  discrete  sum  over  each  of  the  random 
variables  in  the  expression.  The  limits  are  either  the  full  range  of 
the  random  variable  (plus  or  minus  four  times  the  standard  deviation 
is  used  as  an  approximation)  or  are  defined  in  the  expression  being 
evaluated.  Each  range  can  then  be  partitioned  into  a  discrete  number 
(say  100)  of  increments  and  summed  over  all  of  these  increments.  If 
the  increment  has  lower  bound  B1  and  upper  bound  B2 ,  then  the  point 
estimate  (X')  used  in  evaluation  the  expression  is  given  by: 

X'  -  0 . 5*(B2  -  Bl)  (1) 

The  value  for  the  expression  being  evaluated  is  then  multiplied  by 


the  probability  of  each  of  the  point  estimators  occuring,  or: 
Pr {  Bl  <  X'  <  B2  )  -  Pr{X' ) 


(2) 
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This  procedure  is  repeated  in  a  nested  loop  fashion  until  all  of  the 
possible  values  for  each  random  variable  have  been  considered.  The 
summation  expressions  that  are  obtained  using  this  procedure  may 
appear  very  complex,  but  they  are  easily  programmed  as  a  series  of 
nested  loops  with  a  subroutine  to  compute  the  Pr{X')'s  from  equation 
(2).  The  computations  are  straightforward  and  easily  computed  using 
any  personal  computer. 


General  Assembly 

An  overview  of  the  General  assembly  process  is  given  in  Figure 
2.1.  From  Figure  2.1,  it  is  clear  that  the  General  Assembly  process 
doesn't  have  any  assembly  aids  to  limit  variation  in  assembly.  This 
means  that  any  attempt  at  assembly  will  be  successful.  Therefore,  if 
the  completed  assembly  is  within  the  specified  tolerance,  it  will  be 
of  the  required  level  of  quality.  Using  the  nomenclature  given  in 
Figure  2.1,  the  probability  of  the  part  being  within  the  tolerance 
limits  is  given  by  the  intersection  of  the  following  statements: 


Pr {Tx-A< .  5Xy-  ( ex+d*cos  (a- 0)  )<TX+A }  (3) 

Pr{Tx-A<.  5Xy- (ex+d*cos(a+0 )  )<TX+A)  (4) 

Pr {TY-A<.  5YW-  ( cY+d*cos (0-0)  )<TY+A }  (5) 

Pr  {  Ty-  A< .  5YW-  ( €  Y+d*cos  (0+0 ) )  <TY+A }  ( 6 ) 

Statements  (3)  through  (6)  can  be  rewritten  as  follows: 

Pr{Tx-A+ex+d*cos(a-0)<.  5Xw<Tx+A+ex+d*cos  (a- 9 )  }  (7) 

Pr{Tx-A+€x+d*cos(a+0)<.  5Xw<Tx+A+ex+d*cos(a+0 ) }  (8) 
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a  -  cos'x(.5XP/d) 

P  -  cos*1( .5YP/d) 

d  -  [  ( .  5XP)2  +  ( .  5Yp)2  ]0  5 

€x  “  ei  “  0  as  shown  above 

Tx,  Ty  and  A  are  given  constants 

Xp,  X#,  YP,  Yw,  «x,  € y  and  9  are  independent,  normally  distributed 
random  variables  with  known  parameters 


Figure  2.1  :  General  Assembly 
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P  r  { Ty  -  A+  £  Y+d*co  s  (  /?  -  0 )  < .  5  YW<T y+A+  £  Y+d*c  o  s ( £ -  0 )  }  (9) 

Pr(TY-A+£Y+d*cos(£+0)<.  5YH<TY+A+£Y+d*cos  (0+0  )  )  (10) 

Finally,  statements  (7)  through  (10)  can  be  written  as  multiple  sums: 

2  (2  (S  (E(Pr{Tx-A+ex+d*cos(a-  0)<.  5Xw<Tx+A+ex+d*cos(a  -  0)  )  )*Pr  (  0  )  ) 
£xXpYp  0 

*Pr{YP))*Pr{XP))*Pr{ex)  (11) 

2  (E  (2  (2(Pr{Tx-A+£x+d*cos(a+0)<.  5Xw<Tx+A+£x+d*cos(o+5 )  )  )*Pr{  6  ) ) 
£X  XP  Yp  6 

*Pr{Yp))*Pr{XP))*Pr(£X)  (12) 

S  (S  (E  (2(Pr{TY-A+£Y+d*cos(/3-0)<.5Yw<TY+A+€Y+d*cos(/3-0)))*Pr{0)) 
£Y  Xp  YP  e 

*Pr{Yp))*Pr{XP))*Pr{£Y)  (13) 

E  (E  (E  (2(Pr{Ty-A+£Y+d*cos(/J+0)<.5Yw<TY+A+£Y+d*cos(/0+0) )  )*Pr(  9  ) ) 
£Y  XP  Yp  6 

*Pr{YP))*Pr{Xp))*Pr{£Y)  (14) 


These  expressions  are  computed  using  the  numerical  procedure 
described  in  the  previous  section.  Equations  (3)  -  (6)  will  be  used 
several  times  in  the  next  two  sections  and  will  be  abbreviated  as 
T0L1  -  T0L4,  respectively.  In  addition,  it  should  be  noted  that  if 
the  torsional  errors  of  assembly  are  assumed  to  be  insignificant,  the 
results  in  this  and  the  next  two  sections  are  greatly  simplified. 


-in-Hole  Assemble 


An  overview  of  a  Peg- in-Hole  assembly  process  is  given  in 
Figure  2.2.  From  Figure  2.2,  the  probability  that  the  part  can  fit 
into  the  work  is  given  by: 

Pr { Dp<Dh )  (15) 


Likewise,  the  probability  that  the  parts  will  be  assembled  is  given 
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Jf  it 


o  -  cos'l(.5Xp/d) 

0  -  cos'x(.5Yp/d) 
d  -  [  ( .  5XP)2  +  ( •  5YP)2  ]0  5' 

c  -  0.5*(IVDp) 

«,  *  (Y  ■  0  as  shown  above 

Tx,  Ty  and  d  are  given  constants 

XP,  Xy,  YP,  YWt  Dj,,  0p,  «x,  € y  and  8  are  independent,  normally 
distributed  random  variables  with  known  parameters 


Figure  2.2  :  Peg-in-Hole  Assembly 
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by  the  following  two  probability  statements: 


Pr{ - c<  £  x<c ) 

(16) 

Pr{ - C<  6  y<c ) 

(17) 

Finally,  the  probability  that  the  assembled  unit 

is  given  by  the  following: 

is  within  tolerance 

Pr{T0Ll  |  -c<£x<c) 

(18) 

Pr{T0L2  |  -c<£x<c) 

(19) 

Pr { T0L3  |  ~  c< £  y<c ) 

(20) 

Pr{T0L4  |  -  c< £  y<c ) 

(21) 

The  probability  of  manufacturing  an  acceptable  unit  is  equal  to  the 
intersection  of  equations  (15)  through  (21). 

Equations  (15)  through  (17)  are  easily  computed  using  a 
numerical  technique  for  normal  probabilities.  The  last  four 
conditional  probabilities  need  to  be  massaged  into  a  more  usable 
form.  The  probability  that  bottom  corner  of  the  part  in  the 
assembled  unit  will  be  within  tolerance  in  the  X  direction  can  be 
written  as  follows: 


Pr(Tx-A- .  5Xw+d*cos(a-5  )<  -  ex<Tx+A- .  5Xw+d*cos  (a-S)  ,  -c<ex<c  } 

-  (22) 


Pr{  -c<€x<c ) 

The  denominator  of  this  expression  is  equation  (16)  and  is  easily 
computed.  The  numerator  of  this  expression  can  be  rewritten  as 
follows : 


Pr{  -  (Tx+A-  .  5Xw+d*cos(a-5)  )<ex<-  (Tx-A- .  SXw+d^os  (a-  9 ) )  ,  -c<£x<c  ) 


(23) 

Finally,  this  probability  statement  can  be  written  as  the  following 
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multiple  sum: 

2  (E  (E  (  E(Pr{max[ -c ,  -  (Tx+A-  .  SX^dh^cos  (a-$  )  )  ]<ex<min[c  ,  -  (Tx-A- 

C  Xp  Yp.SXy, 

.  5Xw+d*cos(a-  5 ) )  ] )  )*Pr  {  .  SXy)  )*Pr ( Yp)  )*Pr  (XP)  )*Pr  { c } 

(24) 

As  with  the  General  assembly,  this  summation  is  computed  using 
the  numerical  procedure  outlined  in  the  first  section  of  this 
chapter.  The  probabilities  that  the  top  corner  will  be  within 
tolerance  in  the  X  direction  and  both  corners  in  the  Y  direction 
(statements  (19) -(21))  are  computed  in  the  same  manner  outlined 
above.  The  probability  of  producing  a  defective  assembly  is  computed 
in  the  same  manner  used  for  the  General  assembly  process . 

Box- and- Cover  Assembly 

An  overview  of  a  Box  and  Cover  assembly  process  is  presented  in 
Figure  2.3.  Using  the  notation  given  in  Figure  2.3,  the  probability 


that  the  part  can  fit  into  the  work  is  given  by  the  following  two 
probability  statements: 

Pr{XB<Xc)  (25) 

Pr  { YB<YC )  (26) 

Likewise,  the  probability  that  the  parts  will  be  assembled  on  the 
first  attempt  [2]  is  given  by  the  following  probabilities: 

Pr(^<cos'1(.5XB/r)-cos'1((.5XB+dX)/r)  |  dX<dY)  (27) 

Pr(«<cos‘1(.5YB/r)-cos'1((.5YB+dY)/r)  |  dY<dX)  (28) 


where 


dX  -  .5(Xc-Xc)-ex 
dY  -  . 5(YC-YC) - cy 


a  -  cos*1(  .  5XP/d) 

$  -  cos"'( . 5YP/d) 
d  -  [  ( .  5XP)2  +  (  .  5YP)2  ]0-3 
r  -  [  ( .  5Xb)2  +  ( .  5Yb)2  )0  5 
«x  -  «Y  -  0  as  shown  above 
Tx,  Ty  and  d  are  given  constants 

Xpi  Xji,  XB,  Xc,  Yp,  Yh,  Yb,  Yc,  ex,  «p  and  6  are  independent, 
normally  distributed  random  variables  with  known  parameters 

Figure  2.3  :  Box-and-Cover  Assembly 
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Finally,  the  probability  that  the  assembled  unit  will  be  within 
tolerance  is  given  by  the  following  eight  probability  statements: 


Pr {T0L1  |  5<cos*1(.5XB/r)-cos'1((.5XB+dX)/r)  ,dX<dY}  (29) 

Pr{TOLl  |  5<cos'1(  .  5YB/r)  -cos_1(  ( .  5YB+dY)/r) ,  dY<dX)  (30) 

Pr{T0L2  |  5<cos*1(.5XB/r)-cos'1((-5XB+dX)/r)  ,dX<dY)  (31) 

Pr ( T0L2  |  0<cos-1(  .  5YB/r)  -cos'1(  ( .  5YB+dY)/r) ,  dY<dX)  (32) 

Pr { TOL3  |  ^<cos*1(.5XB/r)-cos‘1((.5XB+dX)/r)  ,dX<dY)  (33) 

Pr{TOL3  |  S<cos‘1(.5YB/r)-cos'1((.5YB+dY)/r)  ,dY<dX)  (34) 

Pr{T0L4  |  tf<cos*1( .  5XB/r)  -cos_1(  ( .  5XB+dX)/r)  ,dX<dY)  (35) 

Pr {T0L4  |  S<cos'1(.5YB/r)-cos*1((.5YB+dY)/r)  ,dY<dX)  (36) 


The  probability  of  manufacturing  an  acceptable  unit  is  equal  to  the 
intersection  of  equations  (25)  through  (36). 

The  first  two  equations,  (25)  and  (26),  are  easily  computed 
using  a  numerical  technique  for  normal  probabilities.  Equation  27 
can  be  written  as  follows: 


Pr{0<cos_1(  .  5XB/r)  -cos_1(  (  .  5XB+dX)/r) ,  dX<dY) 

Pr { dX<dY} 


(37) 


The  denominator  of  this  expression  is  easily  obtained  in  the  same  way 

as  equations  (25)  and  (26).  The  numerator  can  be  rewritten  as  the 

following  multiple  sum: 

00 

2  (Z  (Z  (Z  (Z  (  Z  (Pr{J<cos‘1(.5XB/r)-cos'1((.5Xc-€x)/r))) 

Xc  XB  Yc  Yb  «y€X”€xLB 

*Pr(ex))*Pr{cY))*Pr{YB))*Pr{Yc})*Pr(XB))*Pr(Xc)  (38) 

where  exLB  -  .  5(Xc-XB)  - .  5(Yc-YB)+ey 

Equation  (28)  is  almost  identical  to  equation  (27)  and  is  evaluated 


in  a  similar  manner. 
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Statements  (29)  through  (36)  are  all  of  the  same  form  and  the 
evaluation  of  one  sufficiently  illustrates  the  procedure  used  on  all 


eight.  Equation  (29)  can  be  written  as  follows: 

Pr{Tx-A<.  5X#-  (€x+d*cos(a-0  )  )<TX+A,  S<cos_1(  .  5XB/r)  - 
cos_1(  (  .  5XB+dX)/r)  ,  dX<dY) 

-  (39) 


S<cos_1(  .  5XB/r)  -  cos*1  ( ( .  5XB+dX)/r)  ,  dX<dY) 

The  denominator  of  this  expression  is  identical  to  the  numerator  of 

equation  37,  which  has  already  been  computed.  The  numerator  can  be 

rewritten  as  the  following  multiple  sum: 

exUB  SUB 

2  (2  (2  (2  (2  (2  (  2  (  2  (  2  (Pr(Tx-A+ex+d*cos(a-S)<.  5X„ 

X(.  XB  Xp  Yc  YB  Yp  €y  exLB  SUB 

<Tx+a’cxtd*coS(a-»))  <40> 

where  cxLB  -  .  5*(Xc-XB)  - .  5*(Yc-YB)+eY 

exUB  -  +® 


SLB  -  -«o 

SUB  -  cos'1(  .  5XB/r)  -cos_1(  (  .  5XB+dX)/r) 

The  multiple  sums  given  in  equations  (38)  and  (40)  ,  unlike  the  other 
multiple  sums  shown  previously,  aren't  readily  computable.  An 
approximate  solution  can  be  obtained  by  using  the  relationships  for 
the  General  assembly  and  truncating  the  distributions  for  ex,  eY  and 
S.  For  example,  the  distribution  for  S  could  be  truncated  between 
±Smax  where  Snax  is  the  maximum  allowable  torsional  error.  Smax  occurs 
when  ex  -  ey  -  0. 


Comparison  to  Simulation 


The  expressions  derived  using  this  analytical  approach  may 
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appear  intimidating  at  first  glance.  However,  they  are  essentially 
the  same  expressions  required  to  accomplish  a  Monte  Carlo  simulation 
of  the  problem.  In  fact,  that  is  precisely  what  Boucher  [2]  does 
with  them.  The  approach  used  in  this  chapter  is  just  a  small 
manipulation  of  information  required  for  a  simulation  analysis,  but 
it  produces  an  analytic  result  and  not  a  confidence  interval  for  the 
parameter  of  interest. 

Since  the  numerical  procedure  used  in  this  method  is  a  discrete 
approximation  of  a  continuous  probability  density  function,  this 
analytic  result  is  not,  strictly  speaking,  exact.  However,  there  are 
many  routinely  used  procedures  that  are  also  inexact  algorithms. 
Commonly  used  procedures  such  as  long  division  and  finding  the  root 
of  an  expression  are  actually  algorithms  that  follow  a  series  of 
computations  until  a  specified  stopping  rule  is  met.  In  light  of 
this,  the  analytic  results  obtained  using  this  method  are  much  more 
useful  than  the  random  results  of  a  simulated  analysis.  The  analytic 
results  for  the  General  and  Peg-in-Hole  assemblies  are  also  readily 
obtained  using  a  personal  computer  and  don't  require  large  amounts  of 
mainframe  computer  time.  Since  most  of  the  work  required  for  using 
this  method  is  accomplished  when  setting  up  a  simulation,  why  find  a 
simulated,  random  answer  when  a  good  analytic  result  is  readily 


available? 
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Optimization  Model 

Since  an  analytical  result  for  the  probability  of  assembling  a 
defective  unit  is  obtainable,  several  useful  applications  can  be 
made.  The  first  application  is  optimally  allocating  limited  capital 
resources  to  reduce  the  variance  of  the  assembly  processes  for  the 
system  under  review.  The  general  objective  of  the  problem  is  to 
minimize  the  yearly  costs  due  to  capital  depreciation  and  production 
of  defective  units. 

For  this  analysis,  it  is  assumed  that  a  discrete  number  of  new 

processes  are  available  to  choose  from  and  each  of  these  has  an 

initial  cost  that  is  inversely  proportional  to  the  process  variance. 

The  new  equipment  is  assumed  to  be  depreciated  over  a  known  time 

period  while  the  equipment  currently  in  use  is  assumed  to  be  fully 

depreciated.  For  this  reason,  the  cost  of  the  current  equipment  is 

assumed  to  be  zero  for  the  model.  The  optimization  is  summarized 

below  as  a  math  programming  problem: 

nl  n2  Hid 

min  Z  Z  (cij/k)*xiJ  +  P*  Z  Xyl  Z  x2J  .  .  .  (  Z  xmj 

i  J  Ji-1  J2-2  j.-l 

*{D(xlj,x2J,  .  .  •  , ) )  .  .  .}  (41) 

subject  to  Z  Z  CijXij  <  investment  amount 

i  J 

Z  Xij-1  for  all  i-l,2,...,m 

j  1  if  process  j  is  selected  at  stage  i 

XU“ 

0  otherwise 

where  i  -  1,2,... tm  assembly  stages 

j  -  1,2,...,!^  processes  for  each  stage  i 
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k  -  depreciation  factor 

P  -  cost  incurred  due  to  a  defective  assembly 
Cy  -  cost  of  implementing  process  j  at  stage  i 
D(x1j,x2j,.  .  .  .x,^)  -  number  of  defectives  produced  for 
a  given  combination  of  processes 

It  is  clear  that  the  objective  of  this  problem  is  nonlinear,  making 
solution  by  the  usual  linear  techniques  impossible.  However,  a 
simple  manipulation  that  transforms  this  type  of  problem  into  a 
straight  linear  (0-1  integer)  program  is  illustrated  in  Chapter  3. 
Solution  is  then  possible  on  a  personal  computer  using  any 
commercially  available  software  package  for  linear  programming 
problems . 


p  Charts 

Another  application  of  the  analytic  results  obtained  using  the 
method  outlined  earlier  in  this  chapter  is  statistically  controlling 
an  assembly  process.  Since  the  analytically  obtained  result  is  a 
proportion,  controlling  the  assembly  process  using  a  p  Chart  is 
appropriate.  The  expected  proportion  of  defectives  found  with  the 
analytic  approach  is  used  as  the  value  for  p0.  The  statistical  basis 
for  this  type  of  control  chart  is  the  following  hypothesis: 

H0  :  P  -  P0 

Hi  :  p  **  p0 

From  Grant  and  Leavenworth  [5],  the  control  limits  for  this  chart  are 
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given  by: 

UCL  -  Po  +  Za/2((p0-(l-p0))/n)°-5  (42) 

LCL  -  p0  -  Za/2(  (Po-  (l-p0)  )/n)°‘5  (43) 

For  a  -  0.0027,  Za/2  is  equal  to  three.  An  example  of  this  procedure 
is  presented  in  Chapter  3. 
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Chapter  3 
NUMERICAL  EXAMPLE 

Problem  Definition 

The  numerical  example  used  in  this  chapter  is  an  assembly 
consisting  of  four  components  that  are  assembled  in  three  stages. 

The  geometry  of  the  problem  is  given  in  Figure  3.1.  Several 
assumptions  are  used  during  the  analysis  of  this  problem.  First,  the 
person  or  machine  performing  the  assembly  operation  initially 
determines  the  midpoint  of  the  part  to  be  joined  to  the  base  stage. 
The  assembler  then  finds  the  alignment  point  on  the  base  stage.  The 
location  of  this  alignment  point  is  obtained  using  one-half  of  the 
nominal  for  the  part  to  be  joined.  The  assembler  then  attempts  to 
join  the  two  parts  along  these  target  lines.  The  assembler  does  this 
task  with  a  repeatability  error,  for  i  -  1,  2,  3.  The  e^s  are 
assumed  to  be  normally  distributed  with  mean  zero  and  variance  a ^  . 
Once  the  parts  make  contact,  they  are  affixed  to  each  other  in  some 
manner  (glue,  weld,  etc.).  Finally,  the  manufacturing  processes  used 
to  manufacture  the  four  components  are  also  assumed  to  be  normally 
distributed  with  known  mean  and  variance.  The  mean  is  assumed  to  be 
equal  to  the  nominal  dimension  of  the  part  and  the  variance  is 
assumed  to  equal  one -sixth  of  the  total  tolerance. 


e  L  -  18  ±  0.212  W  -  24  ±  0.212 

Lt-  14  ±  0.15  Wx-  20  ±  0.15 

La-  12  ±  0.15  W2-  4  ±  0.15 

L3-  6  ±  0.10  W3-  3.8  ±  0.09 

L*-  20.2  ±  0.15  W*-  3.8  ±  0.09 

ei  -  e2  -  <3  “  NID(0.0,0.031622  ) 


Figure  3.1  :  Numerical  Example 
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Determining  the  Expected 
Number  of  Defectives 

The  first  part  of  the  analysis  determine  the  probability  of  a 
defective  assembly  being  produced  using  the  current  processes.  Since 
there  are  many  combinations  that  can  cause  an  assembly  to  be 
defective,  it  is  easiest  to  determine  the  probability  of  making  a 
defective  assembly  using  the  following  relationship: 

Pr{defective  assembly}  -  1  -  Pr{acceptable  assembly) 

The  probability  of  having  an  acceptable  assembly  is  the  intersection 


of  the  following  probability  statements: 

Pr  ( W  -  A<W  i+W2<W+A  }  (45) 

PrtL-A^+Lj/^-e^L^L+A  |  Ax^Lz/2+e1)  (46) 

Pr{L-A<L,+Lj<L+A  |  A^Lj^+ei)  (47) 

Pr  { A2+W3/2  - «  2<W2 }  (48) 

Pr { W-A<A2+W3/2 - e2+W4<W+A  |  A2>W3/2+e2)  (49) 

Pr{W-A<W3+W4<W+A  |  A2<W3/2+£2)  (50) 

Pr{L-A<A1+L2/2-«1+L3+L*/2-e3-A3<L+A  |  A3<L4/2- e3 , A1>L2/2+£1 }  (51) 

Pr{L-A<L2+L3+L4/2-€3-A3<L+A  |  A3<L4/2-c3,A1<L2/2+ei)  (52) 

Pr  { Aj+L^/2  -  e  j+Lg  -  Lx^L4/2+A3+£3 }  ( 53 ) 

The  conditional  probabilities  are  determined  using  numerical 
integration. 


The  probability  statements  that  don't  contain  any  conditional 
statements  are  easily  computed  using  numerical  techniques.  To  use 


numerical  integration,  the  conditional  probabilities  have  to  be 
rewritten  in  a  more  usable  form.  For  example,  equation  (46)  can  be 
rewritten  as  follows: 
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Pr  { L- A<A1+L2/2  -  €  1+L3<L+A ,  A1>L2/2+c  x ) 

-  (54) 


Pr{A1iL2/2+€1) 

The  denominator  of  this  expression  is  easily  computed.  The  numerator 
must  be  manipulated  into  a  single  probability  statement  on  one  random 
variable  as  shown  below: 

Pr{L-A-A1-L3+€1<L2/2<min[L+A-A1-L3+e1,A1-£1]  )  (55) 

This  statement  is  computed  by  numerically  integrating  over  all 
possible  values  for  each  random  variable: 

2  (2  (Pr(L-A-A1-L3+e1<L2/2<min[L+A-AI-L3+£1,A1-£1)  )*Pr(L3)  )*Pr{  £j} 

€i  ^ 

(56) 

The  actual  value  for  the  original  probability  is  given  by  dividing 
this  expression  by  the  original  denominator.  The  listing  of  the 
BASIC  program  used  to  do  the  numerical  integration  in  equation  (56) 
can  be  found  in  the  Appendix.  The  solution  of  all  the  remaining 
conditional  in  this  example  can  be  found  using  the  same  manner  used 
for  equation  (46) .  It  is  helpful  to  combine  random  variables  into 
joint  probability  distributions  wherever  possible  to  reduce  the 
number  of  summations  required  to  obtain  a  solution. 

After  determining  the  probability  that  the  current  process  will 
produce  a  defective  assembly,  the  following  relationship  is  used  to 
determine  the  expected  number  of  defective  assemblies  produced: 

E(#  of  defective  assemblies  produced/given  time)  - 
Pr(defective  assembly )*{ tot  #  produced/time) 
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For  this  example,  the  probability  of  a  defective  assembly  is  0.02541 
and  the  yearly  production  rate  is  10,000  units.  This  yields  254.1 
expected  defective  assemblies  per  year  using  the  current  processes. 

Sensitivity  Analysis 

A  sensitivity  analysis  performed  on  the  assembly  process  in  this 
example  is  accomplished  using  the  design  of  a  complete  26  factorial 
experiment.  The  six  factors  studied  are  the  mean  and  variance  for 
the  repeatability  errors  for  each  of  the  three  assembly  processes  in 
this  example.  The  two  levels  for  each  mean  are  0.01  and  -0.01.  The 
two  levels  for  each  variance  are  0.0004  and  0.0001.  The  results  are 
summarized  in  Table  3.1. 


Table  3.1  :  Sensitivity  Analysis 


effect  -  (Shieh-Slow) /64 

Mi 

-0.00195 

m2 

0.00011 

m3 

-0.00183 

a\ 

\ 

-0.00066 

-0.00009 

-0.00087 

From  Table  3.1,  it  is  clear  that  the  first  and  third  assembly 
processes  have  the  greatest  effect  on  the  probability  of  producing  a 
defective  assembly.  As  may  be  expected,  a  change  in  the  means  has  a 
greater  effect  than  a  change  in  the  variances. 
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Using  the  results  of  this  sensitivity  analysis,  one  important 
conclusion  can  be  inferred.  Priority  for  any  reductions  in  variance 
should  be  given  to  the  first  and  third  assembly  processes.  This 
conclusion  is  confirmed  in  the  next  section  of  this  chapter. 

Optimization  Problem 

The  premise  of  the  optimization  problem  considered  in  this 
section  is  that  the  assembly  considered  is  in  production  with  known 
process  variances.  There  are  a  discrete  number  of  processes  to 
choose  from  to  improve  the  assembly  operations.  The  initial  costs  of 
these  processes  are  inversely  proportional  to  their  variances.  The 
current  and  potential  processes  are  summarized  in  Table  3.2.  It  is 
assumed  that  a  finite  amount  of  capital  is  available  for 
improvements.  It  is  also  assumed  that  each  defective  assembly 
produced  causes  the  company  to  incur  a  cost  of  P  dollars . 

We  can  further  assume  that  a  depreciating  factor  is  applied  to 
the  capital  outlays,  say  one -fifth,  corresponding  to  a  five  year 
straight  line  depreciation  of  the  new  equipment.  Finally,  if  we 
don't  consider  the  time  value  of  money,  the  general  form  of  the 

optimization  problem  is: 

3  4  4  4  4 

min  E  S  (clJ/5)*xiJ  +  P*  E  x:J  {  E  x2j  {  E  x3J*{D(xXJ  ,x2j  ,x3  ) ))) 

i-1  J-l  Ji-1  j2-l  j3-l 

4 

subject  to  E  xXJ  -  1 

Ji-1 

4 

E  x2J  -  1 

j2-l 


Table  3.2  :  Current  and  Potential  Processes 


current : 
candidate  1 
candidate  2 
candidate  3 


current : 
candidate  1 
candidate  2 
candidate  3 


current : 
candidate  1 
candidate  2 
candidate  3 


Process  1 

cio  “  $0 
cn  -  $2000 
c^2  ™  $4000 
c13  -  $6000 


Process  2 

C20  “  $0 
c21  -  $2000 
c22  *  $4000 
c23  “  $6000 


Process  3 

c30  “  $® 

c31  -  $2000 
c32  “  $4000 
c33  ”  $6000 


a  -  0.03162 
a  -  0.02236 
a  -  0.01581 
o  -  0.01291 


a  -  0.03162 
a  -  0.02236 
a  -  0.01581 
a  -  0.01291 


a  -  0.03162 
a  -  0.02236 
a  -  0.01581 
a  -  0.01291 


4 

2  x3J 

Js-l 
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1 

3  4 

2  2  CijXij  <  investment  amount 

i-1  j-1 

1  if  process  j  is  selected  at  stage  i 

xiJ  “ 

0  otherwise 

This  problem  can  be  made  into  a  straight  integer  programming  problem 
by  introducing  some  constraints  to  ensure  that  one  and  only  one 
process  is  chosen  as  follows: 

P*xn*x21*x3i*  ( D  ( xn ,  x21 ,  x31 ) 
is  rewritten  as: 

P*xd*{D (process  combination  d) } 

xd*xu 

x«j<x21 

x<^x31 

xd5:x11+x21+x31-2 

The  only  drawback  to  this  problem  is  that  it  requires  total 
enumeration  of  all  of  the  possible  expected  number  of  defectives. 
Although  this  sounds  monumental,  it  is  fairly  easily  done  since  the 
basic  probability  statements  are  functions  of  one  or  two  of  the 
assembly  processes. 

A  parametric  analysis  of  this  problem  is  performed  for  several 
levels  of  "defective  costs"  and  investment  capital.  The  results  are 
summarized  in  Table  3.3.  The  results  of  this  analysis  confirm  the 
sensitivity  analysis  performed  earlier.  It  is  optimal  to  reduce  the 


OOOUM* 


27 


Table  3.3  :  Parametric  Analysis 
Cost  per  Defective  Assembly  (P) 
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where  each  cell  represents 
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variance  of  the  first  and  third  assembly  processes  before  the  second 
process.  It  is  also  clear  from  this  analysis  that  it  may  not  be 
advisable  to  expend  all  of  the  capital  available  for  process 
improvement,  depending  on  the  cost  incurred  per  defective  assembly. 
This  analysis  could  also  be  used  to  justify  a  request  for  increased 
availability  of  capital  improvement  funds. 

p  Charts 

A  sample  p  Chart  is  presented  in  this  section  to  illustrate  the 
implementation  suggested  at  the  end  of  Chapter  2.  The  premise  of 
this  illustration  is  that  a  daily  sample  of  50  assemblies  is 
inspected  and  the  proportion  of  defectives  is  recorded  on  the  control 
chart.  Using  the  expected  number  of  defectives  found  earlier  in  this 
chapter  (-  0.0254)  as  pD  and  an  a  of  0.0027,  the  control  limits  are 
found  using  equations  (42)  and  (43).  These  control  limits  are  given 
below: 

UCL  -  0.0922 

LCL  —0.0414  or  0.0 

A  simulated  manufacturing  process  is  run  for  20  days.  For  the  first 
ten  days,  all  of  the  fabrication  and  assembly  processes  are  kept  in 
control.  On  the  eleventh  day,  a  0.06324"  shift  in  the  mean  of  the 
first  assembly  process  (i.e.,  e:  was  now  distributed  N(0 . 06324 , 0 . 001) 
instead  of  N(0 . 0 , 0 . 001) )  is  introduced.  This  shift  remains  in  effect 
through  the  20th  day.  The  p  Chart  for  this  time  period  is  shown  in 
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Figure  3.2.  From  Figure  3.2,  this  simulated  process  is  detected  as 
being  out  of  control  on  the  15th  day.  Once  again,  this  result  is  not 
unique,  since  each  simulation  run  is  different.  This  chart  is 
presented  solely  as  an  illustration  of  the  method  developed  in  this 
paper. 
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Chapter  4 

CONCLUSIONS  AND  RECOMMENDATIONS 


A  method  for  analytically  determining  the  acceptance 
characteristics  of  a  complex  assembly  process  was  developed.  It  was 
shown  that  the  basic  relationships  needed  for  this  analysis  were  also 
needed  for  the  development  of  a  simulation  model,  so  this  method  did 
not  require  much  more  effort  from  the  analyst  than  that  required  for 
a  simulated  result.  The  results  of  this  analytic  method  were  shown 
to  be  very  useful  for  optimally  allocating  limited  capital 
improvements  and  for  assembly  process  control.  An  example  problem 
was  presented  that  illustrated  this  analytic  methodology  as  well  as 
the  two  applications  that  were  developed. 

There  are  many  areas  for  further  research  in  the  area  of 
assembly  process  allocation  and  control.  One  possible  area  for 
further  research  is  the  optimization  problem.  There  may  be  some  way 
to  formulate  the  problem  in  a  continuous  form  that  is  solvable  by 
some  nonlinear  method.  The  problem  can  also  be  expanded  to  include 
the  fabrication  processes  as  well  as  the  assembly  processes,  thus 
providing  a  more  global  optimization  problem  for  a  given  system.  The 
objective  function  may  include  many  more  factors.  For  example,  if  a 
candidate  assembly  process  involves  a  piece  of  special  tooling  that 
not  only  reduces  assembly  variance  but  is  also  more  efficient  for  the 
person  doing  the  assembly,  the  reduced  labor  costs  may  be  added  to 
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the  objective.  Finally,  the  structure  of  the  modified  integer 
programming  problem  is  similar  to  a  network  type  of  problem.  If  the 
constraint  matrix  can  be  manipulated  to  make  the  problem  into  a 
network  type  of  problem,  a  faster  solution  may  be  possible. 

The  Box-and-Cover  problem  also  has  many  possibilities  for 
further  research.  If  some  of  the  random  variables  in  this  problem 
can  be  combined  to  form  joint  distributions,  the  numerical 
integrations  in  equations  (38)  and  (40)  may  become  feasible.  This 
eliminates  the  need  for  an  approximate  formulation  of  the  problem. 
Also,  the  possible  forms  for  the  approximate  solution  of  thi^  problem 
(i.e.,  using  truncated  distributions  in  the  General  assembly  case) 
can  be  further  researched  to  provide  more  accurate  results . 

The  final  area  for  possible  research  is  to  expand  this  method  to 
electronic  and  electromechanical  types  of  assemblies.  Since  purely 
mechanical  assemblies  are  only  one  small  segment  of  industry  today,  a 
method  that  can  be  used  for  different  types  of  assembly  processes  can 
be  very  useful.  For  example,  the  expected  number  of  defectives  for 
an  electronic  assembly  might  possibly  be  used  to  predict  the  expected 
yield  of  the  assembly  through  an  acceptance  test  procedure.  This 
expected  yield  can  then  be  used  to  statistically  control  the  assembly 
process.  This  eliminates  the  need  for  arbitrary  target  levels  for 
acceptance  yield  that  are  commonly  used,  even  though  they  have  no 
statistical  basis.  Management  of  assembly  processes  is  more 
effective  when  the  foundation  of  that  management  is  statistically 
developed  and  analytically  obtained.  Management  by  arbitrary  goals 
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and  simulated  results  is  always  less  than  optimal. 

One  area  for  additional  study  is  to  develop  a  generalized 
computer  program  that  can  handle  a  variable  number  of  summations . 
Ideally,  the  user  is  able  to  input  the  distribution  parameters  and 
the  limits  of  each  sum  interactively.  The  program  should  be  written 
in  a  compilable  language  to  help  speed  the  numerous  calculations.  In 
addition,  the  program  should  minimize  the  number  of  times  that  the 
normal  probability  subroutine  is  called  from  within  a  summation. 

This  can  be  accomplished  through  a  clever  use  of  arrays.  The  program 
in  the  Appendix  is  written  in  BASIC  since  every  computer  has  some 
form  of  BASIC  available  with  the  initial  software. 
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Appendix 

COMPUTER  PROGRAM  USED  TO  PERFORM 
NUMERICAL  ANALYSIS 


20  REM  This  part  of  the  program  is  used  to  compute  std  normal  probs 

21  REM 

25  DIM  Y(6) ,  C(6) 

35  DEFDBL  A-H 
40  DEFDBL  O-Z 

44  REM  Function  to  implement  change  of  variable 

45  DEF  FND  (YT)  -  BO  +  B1  *  YT 

49  REM  std  normal  function 

50  DEF  FNF  (W)  -  1!  /  (EXP((W  A  2)  /  2)  *  SQR(2  *  3.1415927#)) 

55  A  -  0 

60  REM  Input  of  weighting  factors  for  Gauss -Legendre  Formula 
65  FOR  I  -  1  TO  3 
70  READ  C(I) 

75  C(7  -  I)  -  C(I) 

80  NEXT  I 

85  DATA  .171324492, .360761573, .467913935 
90  REM  Input  of  function  arguments 
95  FOR  I  -  1  TO  3 
100  READ  Y(7  -  I) 

105  Y(I)  -  -Y(7  -  I) 

110  NEXT  I 

115  DATA  .932469514, .661209386, .238619186 

116  REM  all  of  the  above  statements  are  used  in  the  subroutine 

117  REM  now  start  the  numerical  integrations 

118  REM 

120  REM  input  parameters  of  given  distributions 

125  INPUT  "mul  ;  MU1 

130  INPUT  "sigmal  ;  SIGMAl 

135  INPUT  "mu2  ;  MU2 

140  INPUT  "sigma2  ;  SIGMA2 

145  INPUT  "mu3  ;  MU3 

150  INPUT  "sigma3  SIGMA3 

160  REM  initialise  suml  and  sum2 

165  SUM1  -  0 

170  SUM2  -  0 

175  L  -  18 

180  W  -  24 

181  A1  -  6 

182  A2  -  1.9 

183  A3  -  1.9 


185  DELTA  -  .212 
200  REM 

210  REM  compute  limits  on  first  summation 

211  REM 

215  UPLIM1  -  4  *  SIGMA1  +  MU1 
220  L0LIM1  -  -4  *  SIGMAl  +  MUl 
225  INCR1!  -  (UPLIM1  -  LOLIM1)  /  50 

230  REM 

231  REM  compute  limits  on  second  summation 

232  REM 

235  UPLIM2  -  4  *  SIGMA2  +  MU2 
240  LOLIM2  -  -4  *  SIGMA2  +  MU 2 
245  INCR2 !  -  (UPLIM2  -  L0LIM2)  /  50 

250  REM 

251  REM  initiate  first  summation 

252  REM 

255  BOT1  -  LOLIM1  -  INCR1 

260  FOR  M  -  1  TO  50 

270  BOT1  -  BOT1  +  INCR1 

280  TOPI  -  BOT1  +  INCR1 

290  PT1  -  (BOT1  +  TOPI)  /  2 

300  B  -  (BOT1  -  MUl)  /  SIGMAl 

310  GOSUB  2000 

320  PR1BOT  -  CUM 

330  B  -  (TOPI  -  MUl)  /  SIGMAl 

340  GOSUB  2000 

350  PR1TOP  -  CUM 

360  PR1  -  PR1TOP  -  PR1BOT 

370  IF  PR1  <  0  THEN  PR1  -  0 

375  REM 

376  REM  initiate  second  summation 

377  REM 

380  B0T2  -  L0LIM2  -  INCR2 

390  FOR  N  -  1  TO  50 

400  BOT2  -  B0T2  +  INCR2 

410  T0P2  -  BOT2  +  INCR2 

420  PT2  -  (BOT2  +  TOP2)  /  2 

430  B  -  (B0T2  -  MU2)  /  SIGMA2 

440  GOSUB  2000 

450  PR2B0T  -  CUM 

460  B  -  (T0P2  -  MU2)  /  SIGMA2 

470  GOSUB  2000 

480  PR2T0P  -  CUM 

490  PR2  -  PR2T0P  -  PR2B0T 

500  IF  PR2  <  0  THEN  PR2  -  0 

505  REM 

506  REM  calculate  probability  for  inner  expression 

507  REM 

510  REM  first  calculate  the  lower  bound 

511  REM 
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515  B0T3A  -  L  -  DELTA  -  Al  -  PT2  +  PT1 

520  BOT3B  -  BOT3A 

524  BOX 3  -  BOT3A 

528  IF  (BOT3A  -  BOT3B)  <  0  THEN  BOT3  -  BOT3B 

530  REM 

531  REM  then  calculate  the  upper  bound 

532  REM 

535  T0P3A  -  L  +  DELTA  -  Al  -  PT2  +  PT1 

540  TOP3B  -  Al  -  PTl 

545  T0P3  -  TOP 3 A 

547  IF  (T0P3B  -  T0P3A)  <  0  THEN  TOP3  -  TOP3B 

550  B  -  (BOT3  -  MU3)  /  SIGMA3 

560  GOSUB  2000 

570  PR3B0T  -  CUM 

580  B  -  (TOP3  -  MU3)  /  SIGMA3 

590  GOSUB  2000 

600  PR3T0P  -  CUM 

605  PR3  -  PR3T0P  -  PR3B0T 

610  IF  PR3  <  0  THEN  PR3  -  0 

620  SUM2  -  SUM2  +  PR3  *  PR2 

630  NEXT  N 

640  SUM1  -  SUM1  +  SUM2  *  PR1 

641  SUM2  -  0 

645  PRINT  "m  — " ;  M 
650  NEXT  M 

660  PRINT  "suml  ;  SUM1 
670  STOP 

2000  REM 

2001  REM  Subroutine  calculates  the  area  under  a  std  normal  curve 

2002  REM 

2005  IF  B  <-  -4  THEN  CUM  -  0 

2010  IF  B  >-  4  THEN  CUM  -  1 

2015  IF  B  <-  -4  OR  B  >  4  THEN  GOTO  2095 

2020  A  -  0 

2025  SUM  -  0 

2030  IJ  -  0 

2035  IF  B  >  0  THEN  GOTO  2055 

2040  IF  B  -  0  THEN  GOTO  2090 

2045  B  -  -B 

2050  IJ  -  1 

2055  BO  -  (A  +  B)  /  2 

2060  B1  -  (B  -  A)  /  2 

2065  FOR  I  -  1  TO  6 

2070  SUM  -  SUM  +  C(I)  *  FNF(FND(Y(I) ) ) 

2075  NEXT  I 

2080  SUM  -  B1  *  SUM 

2085  IF  IJ  -  1  THEN  SUM  -  -SUM 

2090  CUM  -  .5  +  SUM 

2095  RETURN 


